
rm(list = ls()) #清空内存函数
library(RMySQL)
con <- dbConnect(MySQL(),host="192.168.3.139",
                 dbname="niek_160906",user="niek",password="FbWf8AKC")  
dbSendQuery(con,'set names utf8')

##################################################################################
#设定指数统计周期（xx月~xx月）
date_begin<-24133  #最早为2009年1月
date_end<-24240

#目标产品类型
pty<-11

#基期及城市设定
year_base<-2015
date_base_begin<-24181
date_base_end<-24192
date_base<-24133  #基期设定为2011年1月
city_code<-'qd'   #城市设定为青岛
pty_base<-11      #基本产品类型设定为住宅
pts_base<-10000   #基点设定为10000
################################################################################################################

#数据抽取
source('R/load_spatial_data.R')

#ha_price与ha_phase数据处理
library(dplyr)
price_sale<-subset(tabln.vec$ha_price,proptype==pty_base)[,c('ymid','ha_code','saleprice','salebldgarea')]%>%na.omit
ppi_sale<-merge(price_sale,ha_info.sp,by='ha_code')

data_base<-subset(ppi_sale,ymid>=date_base_begin&ymid<=date_base_end&city_code==city_code)

myvar_city<-names(data_base) %in% c('ymid','city_code','ha_code','proptype','salecount','name','ha_cl_code','ha_cl_name','dist_name','volume_rate','greening_rate')
if(nrow(data_base)>=30){data_base2<-data_base[!myvar_city]%>%
    group_by(dist_code)%>%
    summarise_each(funs(mean))%>%as.data.frame
}else{stop("该基期样本量不足30，需重新设定")}

#设定row.names
data_base3<-as.data.frame(data_base2,row.names = as.character(data_base2$dist_code))%>%
    select(-one_of('dist_code'))

dist_nm<-as.character(data_base2$dist_code)
dist_len<-length(dist_nm)

#基期时点基点为10000的估价计算
library(spdep)
data_base4<-subset(ppi_sale,ymid==date_base&city_code==city_code)
myvar_city2<-names(data_base) %in% c('ha_code','ymid','city_code','proptype','salecount','name','ha_cl_code','ha_cl_name','dist_code','dist_name','volume_rate','greening_rate')
data_base5<-data_base4[!myvar_city2]

#计算各行政区的样本量
dist_num<-table(data_base4$dist_code)%>%as.data.frame
names(dist_num)<-c('DIST','Freq')

#建立权重矩阵
sp_kn<-cbind(data_base5$x,data_base5$y)%>%
    knearneigh(k = 5)%>%
    knn2nb(row.names=rownames(data_base5))

sp_w<-nb2listw(sp_kn,style="W", zero.policy=TRUE)

# model building
data_base6<-data_base5%>%select(-one_of('x','y'))
fit_base<-lagsarlm(log(saleprice)~., data= data_base6, listw=sp_w,
                   zero.policy=T,type = 'Durbin',tol.solve = 1.0e-16)

#predict
data_new<-rbind(data_base3,data_base5)

sp_kn_new<-cbind(data_new$x,data_new$y)%>%
    knearneigh(k = 5)%>%
    knn2nb(row.names=rownames(data_new))

sp_w_new<-nb2listw(sp_kn_new,style="W", zero.policy=TRUE)

data_new2<-data_new%>%select(-one_of('x','y'))
price.pre<-predict(fit_base,newdata=data_new2,listw=sp_w_new,pred.type = 'TS',legacy.mixed = T)%>%
    head(dist_len)%>%exp


pre_dist<-as.data.frame(cbind(dt<-rep(date_base,dist_len),dist_code<-dist_nm,pre<-price.pre))
names(pre_dist)<-c('DT','DIST','PRE')

pre_dist_final<-merge(dist_num,pre_dist,all.x=TRUE)

price.pre_base<-sum(pre_dist_final$Freq*as.numeric(as.character(pre_dist_final$PRE)))/sum(pre_dist_final$Freq)

##data_base3为基期的标准建筑各特征值
##pre_dist_final为各行政区的预测价格
##price.pre_base为该城市在基期的预测价格
